@rdkoster, @biljanaorescanin:

I started looking into Biljana's latest 1-day AMIP simulations using the "develop" branch ("stock") and this PR's branch ("exp").
I didn't quite see what I was expecting before I ran out of time.

However, it looks like the apportioning of the EVACC term in GEOS_CatchGridComp.F90 when running in coupled land-atm (GCM) mode has two errors. The first error was that the EVACC correction was only applied for SUMEVP>0, which should already fixed in the branch. The second error seems to be a confusion of units. For example, in the following equation:

EVPICE = EVPICE - EVACC*EVPICE/SUMEV

EVPICE and SUMEV are in W/m2 and EVACC is in kg/m2/s. It looks like a units conversion with L_s (latent heat of sublimation, MAPL_ALHS) is missing. Similarly, the next three equations are missing a units conversion with L_v (latent heat of vaporization, MAPL_ALHL). In the current code, the correction is therefore essentially ignored (besides being only applied for evaporation but not for dewfall) and the EVP[*] output becomes inconsistent with the EVAP and LHLAND variables.

The error might have been prevented if the code included clear documentation of the units of the inputs and outputs of catchment(). This information is sorely needed to make the code understandable to future investigators.

If the above makes any sense, please modify the branch accordingly and re-run the 1-day AMIP experiments.

@biljanaorescanin: To be safe, it would be best to switch to a restart date in mid-February. We tend to have the biggest problems with the snow terms, and the snow-covered land area is larger in Feb than in Apr, so any issues with snow will be easier to detect.

Thanks, @biljanaorescanin, for making the proposed changes and running another day in Feb.

For this PR's branch at d95a4a1, the following energy balance:

ECHANGE = SWLAND + LWLAND - SHLAND 
          - LHLAND 
          - Lf*PRECSNO - SPLAND - SPSNOW - SPLH

closes within +/- 0.2 W/m2 for all grid cells, with a std-dev across all grid cells of ~0.02 W/m2.
Moreover, we now have

EVLAND = (EVPINTR + EVPSOIL + EVPTRNS)/Lv + EVPSBLN/Ls

to within 1e-3 mm/d. [CORRECTED ON 12 JULY 2024, ORIGINAL NUMBER OF 1e-5 mm/d WAS A TYPO.]
Alas, we still only have

LHLAND = EVPINTR + EVPSBLN + EVPSOIL + EVPTRNS

to within +/- 2 W/m2. The mismatch seems to occur only when FRSNOW>0, but I'm not sure. This requires more investigation.

For reference, in "develop" the following energy balance:

ECHANGE = SWLAND + LWLAND - SHLAND 
          - EVPINTR - EVPSBLN - EVPSOIL - EVPTRNS
          - Lf*PRECSNO - SPLAND - SPSNOW 

is equivalent to the branch's balance above (i.e., closes to within +/- 0.2 W/m2 for all grid cells, with a std-dev across all grid cells of ~0.02 W/m2).

However, for "develop":

EVLAND = (EVPINTR + EVPSOIL + EVPTRNS)/Lv + EVPSBLN/Ls

does not match (deviates up to 10 mm/d) and

LHLAND = EVPINTR + EVPSBLN + EVPSOIL + EVPTRNS

does not match (deviates up to 400 W/m2).

cc: @rdkoster

0-diff tests for GEOSldas successfully completed by @biljanaorescanin on 8 July 2024

@gmao-rreichle last commit is @rdkoster suggestion how to change units. I will send run output once done in email.

@gmao-rreichle last commit is @rdkoster suggestion how to change units. I will send run output once done in email.

@rdkoster, @biljanaorescanin: Thanks again for implementing the change. At first glance it does make more sense to do the corrections in this way. (As an aside, in the new corrections, we should be able to just use EVLAND instead of computing SUMEV; likewise, we could probably have used LHOUT instead of computing SUMEV in the old corrections).

Alas, the change does not solve the remaining issues. To recap, the issues are:

1. LHLAND=EVPSUM is true only to within +/- 2 W/m2 at a subset of tiles.

2. When using LHLAND in the energy balance, the balance only closes to within +/- 0.2 W/m2 at a subset of tiles. (The balance is worse -- to within +/- 2 W/m2 -- if we use EVPSUM instead of LHLAND.)

The latest commit only changes EVPSUM and not LHLAND, the energy balance using LHLAND is unchanged. When using the revised EVPSUM in the energy balance, the balance is no better than with old EVPSUM (to within +/-2 W/m2). Also, with the revised EVPSUM, the equality EVPSUM=LHLAND is still true only to within +/- 2 W/m2.

It's probably still true that the lack of closure only happens where FRSNOW>0, although there are tiles with FRSNOW>0 that do have closure. But I ran out of time to look into this further. I might also have messed up my plots in the rush to take a quick look. We'll have to revisit next week.

@zhaobin74 : I changed three long names of CICE longwave exports in commit feb6d1d
These changes seem to have been missed in #764, but I'm not sure. Can you please double-check and let me know?
Note that the present PR is on top of develop. I have another PR in the works that is on top of R21C. The latter is urgent. I would greatly appreciate it if you could take a look soon. Thanks!

@zhaobin74 : I changed three long names of CICE longwave exports in commit feb6d1d These changes seem to have been missed in #764, but I'm not sure. Can you please double-check and let me know? Note that the present PR is on top of develop. I have another PR in the works that is on top of R21C. The latter is urgent. I would greatly appreciate it if you could take a look soon. Thanks!

@gmao-rreichle, just back from a long leave. I missed changes to these variables in #764. Thanks for adding them.